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1. Introduction 

Low-density parity-check codes (LDPC) have attracted significant interest in recent 
years due to their simphcity and exceptionally high performance I.j. Their simplicity 
and inherent randomness make them amenable to analysis using established methods in 
the area of statistical physics. These have been employed in a number of papers [21-10] 
to gain insight into the properties of LDPC codes and to evaluate their performance. 

These studies include the evaluation of critical noise levels for given codes [2], 
an exact calculation of weight and magnetisation enumerators jlj, the performance of 
irregular codes 0, properties of codes in real- valued channels and the derivation of 
bounds for the reliability exponent P, to name but a few. These studies also represent 
the interdisciplinary nature of this research area and illustrate the successful interaction 
between researchers in the two disciplines. 

The evaluation of error exponents has been a long-standing problem in information 
theory pUllllj . Efforts to obtain exact expressions and/or bounds to the error exponent 
resulted in partial success; although tight bounds have been derived in the case of 
random codes and LDPC with infinite connectivity |Tn|, only limited results have been 
obtained for sparely connected codes. Main stream techniques to tackle the problem 
include sphere-packing and union-bound arguments EJj- Below a certain code- 
rate value, the estimated bounds also become loose and require using the 'expurgated 
exponent' techniques PO] for obtaining a tighter bound. 
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In this paper, we employ methods of statistical physics to evaluate directly the 
average error exponent and typical reliability exponent in Gallager and MN LDPC 
codes. The average error exponent is obtained by carrying out averages over the 
ensemble of randomly generated LDPC codes of given rate and connectivity; while 
the reliability exponent is obtained by selecting the best codes in that ensemble. 
Averages result in the emergence of macroscopic variables, representative of the ensemble 
properties, that can be obtained numerically and used to calculate the average error 
exponent (in the current calculation we assume that short loops, which contribute 
polynomially to the block error probability in LDPC codes ITS], have been removed). 
Average error exponent solutions have been obtained for both finite and infinite 
connectivity vector ensembles, while reliability exponent solutions have been obtained 
only in the case of infinite connectivity. 

As a reference point to test our theory, we use known results obtained in the 
information theory literature for solvable limits (e.g. codes of infinite connectivity), 
and find that our method reproduces them exactly. Perhaps not surprisingly, we also 
find that at fixed noise level and code rate, the reliability exponent for codes of finite 
connectivity is always upper-bounded by that of the infinite-connectivity case. 

Before we proceed, the distinction between the typical bounds found previously 
using methods of statistical physics 0, and the current calculation should be clarified. 
In the former, one employs methods of statistical physics to calculate the typical value 
of a bound based on inequalities introduced by Gallager; while in the current calculation, 
a direct estimation of the average error exponent, rather than a bound, is sought. An 
additional advantage of the current approach is that it can be extended to provide 
reliability exponent values for LDPC codes by restricted averages over codes of high 
performance. 

The paper is organised as follows: In section |2l we introduce the general coding 
framework and the technique used. In sections 01 and ID we present an outline of the 
derivation and the solutions obtained in both finite and infinite connectivity cases 
respectively. In Sectional we compare the error exponent results obtained for MN codes 
to those of Gallager codes in both finite and infinite connectivity cases. Discussion and 
conclusions are presented in section IHl 

2. Definitions 

A regular {k,j) Gallager error-correcting code is defined by the binary (A^ — K) x 
(parity check) matrix A = [C1IC2], which is known to both sender and receiver. The 
(N — K) X [N — K) matrix C2 is taken to be invertible. The number of non-zero elements 
in each row of A is given by k, while the number of non-zero elements per column is 
given by j = k{N - K)/N. 

Gallager's encoding scheme consists of generating a codeword t G {0, 1}^ from 
an information (message) vector s G {0, 1}^ (with A^ > K) via the linear operation 
t = G^s (mod 2) where G is the generator matrix defined by G = [/[Cg^^Ci] (mod 2). 
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The code rate is then given hj R = K/N = 1 — j/k, and measures the information 
redundancy of the transmitted vector. 

Upon transmission of the codeword t via a noisy channel, taken here to be a BSC, 
the vector r = t-\-n^ (mod 2) is received, where G {0, 1}^ is the true channel noise. 
The statistics of the BSC is fully determined by the flip rate p G [0, 1]: 

P(n°) = (1 -p)(5„o,o+P^nO,i (1) 

Decoding is carried out by multiplying r by A to produce the syndrome vector 
z = Ar = y4n°, since AG^ = by construction. In order to reconstruct the original 
message s, one has to obtain an estimate n for the true noise n^. First we select 
the parity check set of A and nO, i.e. all n that satisfy the parity check equations: 
Jpc(v4,n°) = {n \ An = An°}. Since all operations are performed in modulo 2 
arithmetic, Tpc(v4, n°) typically contains exp[A^i? ln(2)] candidates for the true noise 
vector n°. 

It was shown (see e.g. [3 El El for technical details) that this problem can be 
cast into a statistical mechanics formulation, by replacing the field ({0, 1}, +mod(2)) 
by ({1,— l},x), and by adapting the parity checks correspondingly. From the parity 
check matrix A we construct the binary tensor A = 1 <ii <i2 ■ ■ ■ <ik^ N}, 

where Ai^i^.-.i,^) = 1 ii A has a row in which the elements {ic-, c = 1, . . . , /c} are all 1 (i.e. 
when the bits {ii ■ ■ - ik) are involved in the same parity check), and otherwise. The 
fact that each bit ii = 1, . . . , is involved in exactly j parity checks is then expressed 
by X]i2<-<ifc -^(n-'^fc} — h V zi = 1,...,A^ and the parity check equations become 

Decoding now consists in selecting an n from Tpc(.4., n°), on the basis of its 
noise statistics, which are fully described by its magnetisation m{n) = l/NY^^Ui 
(corresponding to the weight in the information theory literature). Note that the number 
of flipped bits in a candidate noise vector n is given by A^(l — m(n))/2. Therefore, we 
introduce a Hamiltonian or cost function for each noise candidate that is negatively 
proportional to its magnetisation: 

N 

H{n) = -F^ni = -FNm{n) (2) 

i=l 

where we take F = ^log^^, such that up to normalisation exp(— /7(n)) yields the 
correct prior for candidate noise vectors generated by the BSC Then, a vector n 
from Ipc{A, n°) with the highest magnetisation (lowest weight) is selected as a solution; 
this corresponds to Maximum A Posteriori (MAP) decoding. 

We are now interested in the probability that other candidate noise vectors are 
selected from the parity check set Xpc{A,n^), other than the correct (i.e. true) noise 
vector n°, for any given combination {n°, A\\ this is termed the block error probability. 
In order to calculate this probability, we introduce an indicator function: 

A(nO,^) = hm hm [Z^(n°,^;A) Z^,\n\ A, (3) 
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where 

Z,in',A;P,)= Yl Z,{n\ A, (3,) = e-^-^^^'l {A) 

neipcino ,A)\no n€ipc{no,A) 

The partition functions Zi{n^, A; f3i) and Z2(n°,A;/?2) differ only in the exclusion of 
from Zi. If the true noise n° has the highest magnetisation of all candidates in the 
parity check set (decoding success), the Boltzmann factor exp[— /5if(n°)] will dominate 
the sum over states in Z2 in the limit f3 — ^ 00, and A(n°,^) = 0. Alternatively, if 
some other vector n ^ -nP has the highest magnetisation of all candidates in the parity 
check set (decoding failure), its Boltzmann factor will dominate both Zi and Z2 and 
A(n°, ^) = 1. Note that the separate temperatures /5i and [32, which are put to be equal 
to (3 in the end, and the powers Ai^2 which are taken to be ±A in the end, have been 
introduced in order to allow us to determine whether obtained solutions are physical or 
not. The power A > have been introduced to restrict the indicator function results to 
0/1. In principle, this can be done by taking the limit A ^ 0; however, in section IHl we 
show that finite < A < 1 values will be used due to various constraints. 

To derive the average error exponent, we take the logarithm of the indicator function 
average with respect to all possible realisations of true noise vectors n*^, and the ensemble 
of regular (fc, j) codes A: 



where 



and 



ir/AW nn = l '^[Z]i2<--<ife ■^(jl---«fc> j] /(■^) 

{t{A))A= ^ .r^ • V) 

To obtain an expression for the reliability exponent one carries out a similar 
calculation with one main difference: prior to averaging the indicator function over 
the ensemble of regular {k,j) codes A, one takes the averaged expression with respect 
to realisations of true noise vectors n° to a power r which favours code constructions 
with a low average error probability (i.e., r < 1). The logarithm of the expression 
averaged over the ensemble of codes A is then divided by r to remove the exponent. 
The expression calculated is: 

Qr = i^^^,^og{[{A{n\A))^.,Y)^ (8) 

Since there are only discrete degrees of freedom, physically meaningful solutions must 
have a non-negative entropy, requiring the disorder-averaged entropies of the two 
partition functions Q to be non-negative. Note that due to the order of taking the 
logarithm vs the various averages, expressions © and (jHl) are not equivalent to a 
(quenched) disorder-averaged free energy. Using general principles one can show that for 



Average and reliability error exponents in low-density parity-check codes 



5 



general values of Pi^2 and Ai^2, the disordered-averaged entropies (with averages taken 
over the joint distribution of code-constructions {A}, true- and candidate-noise {n^, n} 
as suggested by ^ and (jH))) are given, for both calculations © and (jH)), by 



which have to be positive. 



X = 1,2 



(9) 



3. Average error exponent - general solution 



Using standard statistical physics methods such as in [m], we perform the gauge 
transformation Ui —>■ UiU^, and the averages over true noise (jHl) and code constructions 
([7j). In the case of r 7^ 1, each quantity carries two indices (a replica index and another 
index coming from the power r); however, the two indices factorise unless an explicit, 
more complex, symmetry breaking structure is introduced. Here, we do not assume a 
more complex structure that entangles the two types of indices; we also assume the 
simplest replica symmetric scheme ^Hj to arrive at the following expression for the 
average error exponent (r = 1), and for the reliability exponent (optimised r): 

1 



where 



-Extr^ 



log Ii[tt] -j\og /2[vr,7r] + log /3[7r] 



^ c=l 



1 + nti 



rX2 



{d7r{x,y) d7i{x,y)} 



l+xx\ I 1 + yy 



(10) 



:iii 



(12) 



JJ{rf7r(xc,yc)} 



c=l 



u=±l 



c=l 



1 + UXc 



Ai 



X 



E 



v=±l c=l 



n 



1 + vijc 



Ms 



(13) 



where we have used the short-hand notation df{x, y)=dxdy f{x, y). For r = 1, functional 
extremisation of (fTUj) with respect to the densities tt^x, y) and tt{x, y) results in a closed 
set of equations (reminiscent of 'density evolution' equations P): 



„ fe-i 

T^{.x,y) = / JJ{rf7r(xc,yc)}5 

c=l 



c=l 



6 



k-l 



y-Ylvc 



c=l 



7r{x,y) 



(14) 



(15) 
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where 

= ]J{d^x,,y,)}{Dl^{x-,(3^)D>;:^{y-,P,) ■ )^,, (16) 

c=l 

D^{z; (3) ^ [e^^"" + z,)] ± [e"^^"" - z,)]. (17) 

c=l c=l 

For given (/5i, /?2, ^i, A2) in general, solutions to (fT^ and (fT3j) can only be obtained 
numerically. Inserting these solutions into (fTUj) we then obtain Q(/?i, /?2, ^i, A2), which 
becomes the average error exponent for Ai = — A2 = A > 0, and for /5i = /32 = /5 00. 

We must recall, however, that physically meaningful solutions must satisfy the 
conditions © stating that the entropies related to the full and the restricted partition 
sums are non-negative. 

We restrict ourselves to regions below the thermodynamic transition where the 
average case is dominated by the ferromagnetic solution, such that we can fix the system 
described by Z2 in Q to the ferromagnetic solution. This dominance is guaranteed if 
the following constraint is satisfied 
dQ 
d(3 

It turns out that for given A > 0, the largest value of /5 for which (fTHj) is satisfied is 
given by the simple expression /5 = 1/(1 + A). Hence, in order to maximise (3, we must 
look for the smallest value A* that satisfies the conditions on the non-negativity of the 
entropies Q. Unfortunately, in general this value A=k can only be obtained numerically. 
The value obtained for the average error exponent by this analysis is then given by 
Q(l/(1 + A,), 1/(1 + A,), A„ -A,) from 

In figure Ewe present the obtained average error exponent as a function of the flip 
rate for {k,j) = (4,3), {R = 1/3) and {k,j) = (6,3) {R = 1/2) codes. We observe that 
the error exponent indeed converges to zero, as it should, when the flip rate approaches 
its critical value. 

Notice the similarity between the equations obtained here and in in spite of the 
different starting points. It has been shown in [Hj that the analysis should be refined in 
low rate regions by considering a more complex bound. The refined analysis resulted in 
tight bounds of the error exponent even in the region of low code-rates, similar to those 
obtained using expurgated exponent methods. In the next section we will show that 
the selection of 'best codes' through the optimisation of the power r, in calculating the 
reliability exponent, provides similar results to those obtained in 



< . (1^ 

Ai=— A2=A 



4. An exactly solvable limit: k,j^oo 

Whereas for finite density codes solutions for the average error exponent are obtained 
numerically, in the limit of k,j ^00 (while keeping the rate R = 1 — j/k finite) one 
obtains two types of analytic solutions to equations ()14|) and ()15|1 . which can be verified 
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Figure 1. Average error exponent Q as function of the flip rate p for codes of 
= (4,3) (left picture) and {k,j) = (6,3) (right picture). Lines with markers 
correspond to the finite (fc, j) cases. For comparison we also present (thick solid lines) 
the value of the average error exponent in the case of k,j — > oo with R = 1/4 (left) 
and i? = 1/2 (right) as described in the analysis of section^ Note that the transition 
from type I to type II solution occurs at small p values outside the range of this figure. 



by substitution. Moreover, in this limit one also obtains solutions in the reliability 
exponent calculation (jH)), which are generally difficult to obtain for finite k and j values. 



Average error exponent 

Solutions obtained in the average error exponent calculation take the following form: 
Type I: 

1 



Tr{x,y) = -[6{x-l) + 6{x + l)] 6{y - 1) 
ir{x,y) = ^[Six-l) + S{x + l)] S{y - 1) 



(19) 



Type II: 



Tc{x,y) = 5{y - 1) 



G+{F{l+(32\i)) 5(x-tanh(/3iF)) 



+ G_(F(1+/52A2)) (5(x+tanh(/5iF)) 
Ti{x,y) = 6{y-l)6{x) (20) 

with G±{x) = |[1 ± tanh(x)]. 

Taking Pi = P2 = P and Ai = — A2 = A, the average error exponent as obtained 
from the type I solution is given by 

= log 2 - log cosh F + log cosh(/3FA) + log 2 cosh(F - pFX) . (21) 
k 

We find that the entropies Q are always identically zero, and that the constraint (|TH|l 
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requires that (3 = 1/2, such that A = 1 and 



Qi = log 2 - log[e^ + e-^] + log[e^ + e"^ + 2] 

which is exactly the Bhattacharyya limit [TT] . 

The average error exponent as obtained from the type II solution is given by 



(22) 



Q 



II 



A 



log 2 + log 2 cosh 

rh 



+log[2 cosh(F-/3FA)]-log 2 cosh F(23) 



The condition on the entropy {S2) > is satisfied for all /3 > 0, whereas the condition 
{Si) > is violated below the critical (freezing) temperature obtained from 



- ^ log 2 - f]*F tanh[/5*F] + log 2 cosh[/3*F] = 
k 



(24) 



This negative entropy is an artifact of the assumption about the symmetry between 
replicas, and is easily remedied by considering a 'frozen RSB' ansatz [2 . Using this 
ansatz and taking into account condition (fTHj) . the (frozen) average error exponent 
obtained from the type II solution, is finally given by 



Qf^ = F tanh[/5*F] + ^ log 2 - log 2 cosh F 

rC 



J 



(25) 



What remains is to determine whether the type I or type II solution is physically 
dominant, by using Q as a generating function for calculating the related free energies 
(through its derivative with respect to A). Results for the case of k,j^oo are presented 
in figure El for p = 0.01 and p = 0.05. 



4.2. Reliability exponent 

To obtain the reliability exponent we take equations (jl0|) - (jl3p and optimise with respect 
to r. Deriving a general set of equations similar to p4ll5|) . that can be solved iteratively, 
is difficult in this case. However, in the limit k,j — > 00, we observe that we can restrict 
the possible solutions of 7r{x,y) to two different types: 
Type I: 

7r{x, y) = ^ [6{x - 1) + 6{x + 1)] 6{y - 1) (26) 

Type II: 

n{x,y)= 6ix)6iy-l) (27) 

In this case, knowledge of the solution for 7r(x, y) is sufficient for calculation the 
reliability exponent (fTIUl . Furthermore, the expression obtained from the type II solution 
turns out to be identical to that of the average error exponent ()25|) . 

On the other hand, the reliability exponent obtained from the type I solution is 
somewhat different, and takes the form: 

Qr I = --f log2-log[cosh(F)] + -log[cosh''(F)+cosh"((2/?A-l)F)] .(28) 
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Given the relation (|18|) and /5 = 1/(1 + A), one obtains X = 1, f3 = 1/2, and the 
expression reduces to 

Qr I = --Uog2 - log[cosh(F)] + -log[cosh^'(^) + 1] . (29) 

V K T 

Optimising the expression with respect to r, one obtains a similar expression to the 
expurgated exponent result [10] 

Eex{r,R) = max|ln2coshF- ^ In [(2 cosh F)'' + 1] + ^(1 - i?) In2| , (30) 

which is also identical to the result obtained for the average bound of the reliability 
exponent in [U]. 

The reliability exponent is therefore identical to the average error exponent except 
for very low R values as shown in figure|21for p = 0.01 and p = 0.05 (marked by a dotted 
line in the two cases considered). 




R 



Figure 2. Reliability and average exponents as function of the code rate R for regular 
/c, j — > oo Gallager codes for which analytical expressions can be derived; see l)22|l and 
(|25|1 (dashed: p = 0.01 and solid: p = 0.05). The reliability exponent is identical to 
the average error exponent except for very low R values where it is represented by 
the curved solutions above the linear average exponent results marked by a dashed 
and solid lines, respectively. The transition point is marked by a vertical line. The 
transition between solutions of type I and II is marked by • and the critical transition 
point by a X . 
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5. MN codes 



In this section, we extend our treatment of the average error and rehabihty exponent to 
regular MN codes ^21; a variant of LDPC codes. 

A regular MN code is defined by the binary x (A^ + K) matrix A = [Cs\Cn], 
concatenating two sparse matrices with the N x N matrix C„ assumed invertible. The 
N X K matrix Cg has k non-zero elements per row and j per column while C„ has t 
non-zero elements per row and per column. The code rate is given hjR = K/N = k/j. 
The encoding scheme consists of generating a codeword t° G {0, 1}^ from an (unbiased) 
message vector s° G {0, 1}''^ via t° = {C^^Cs)s^ . Upon sending through the noisy 
channel the vector r = t° + n° is received, where n° is the true channel noise (^. 

Decoding is carried out by multiplying the received vector r by C„ to produce 
the syndrome vector z = CgS^ -\- CnTi^ ■ In order to reconstruct the original message, 
one selects the best estimate (s,^), for the true (s°,n°) from the parity check set 
Xpc = {(s, n) I CgS + C„n = z}, on the basis of the message/noise statistics. Note that 
since we take the message vector s° to be unbiased, the selection will only be based on 
the noise statistics. 

Since most calculation steps are completely analogous (although lengthier) to those 
of Gallager codes, we only state the final general expression for MN codes: 

QriPi, P2, Ai, A2) = -Extr^,^-,p,p I (31) 



— k log 
-tlog 



1 



+ log 



{dn{x,y) dn{x,y)} 

{dp{x,y) dp{x,y)} 
t 

JJ {rf7r(xi, yi)} JJ {dp{uu Vi)} 



XX 



2 

1 + XX 



2 

+ yy 



2 J \ 2 



rX2 



rX2 



1=1 



rXi 



r\2 



+ - log / n{c?vr(xc,yc)} 

c=l 

+iog / n{^^(^''^')} 
1=1 



En 



=± c=l 



1 + axc 



rXi 



,l3iTFn'> 



n 

1=1 

T=± 



1 + TXi 



n 



En 

(T=± c=l 

Ai 



1 + ryi 



r\2 



with the short-hand notation df{x,y)=dxdy f{x,y). As for Gallager codes, for jSi = 
(32 = (3, and Ai = — A2 = A, Qr becomes the average error exponent for r = 1 , while 
for optimised r it becomes the reliability exponent. Furthermore, the conditions (jH)) and 
(fT5|) must always be satisfied. 

Similarly to the case of Gallager codes one can derive a set of functional equations 
(reminiscent of 'density evolution' equations jT]) for 7T,n,p and p. 
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Figure 3. Average error exponent Q as function of the flip rate p obtained for Gallager 
and MN codes Left: results for MN codes of i? = 1/4. Lower line and circles: 
{k,j,t) = (1,4,2). Upper line and diamonds: (2,8,3). Thick upper line: k,j —> oo. 
Right: results for Gallager and MN codes of i? 1/2. Lower line and circles: Gallager, 
{k,j) — (6,3). Upper Hne and diamonds: MN, {k,j,t) — (3,6,3). Upper Une: the 
analytical solution of fc, j oo obtained via (|22I25|I . 



5.1. Average error exponent - finite k,j and t 

The average error exponent can be calculated numerically for finite k,j and t values; 
the average error exponent Q as function of the flip rate p is shown in figure 01 

On the left, we show results for MN codes of fixed rate R = 1/4 with three different 
sets of parameters {k,j,t) = (1,4,2) (circles), (2,8,3) (diamonds) and k,j oo (upper 
line). It is interesting to notice that average exponents for either k > 2 or t > 2 values 
coincide with that of the infinite connectivity case (which can be obtained analytically). 
This complements other interesting properties of MN codes, to do with their critical 
flip rate values, that have been obtained previously, distinguishing them from Gallager 
LDPC codes 

On the right, we see a comparison between average error exponents of Gallager 
and MN codes {R = 1/2). The Gallager code {k,j) = (6,3) average error exponent 
(circles) is significantly below the random code k,j—^oo value (thin upper line) and 
the equivalent MN code {k,j,t) = (3,6,3) result (diamonds). 

5.2. Average and reliability error exponents - k, j,t ^ oo 

The case of k,j,t ^ oo is solvable exactly for all transmission rates, and both average 
and reliability error exponents can be obtained analytically. The solutions obtained 
as well as the average and reliability error exponents calculated are identical to those 
of Gallager LDPC codes. Retrospectively, this is not surprising as both codes become 
random codes in this limit. 
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6. Discussion 

In this paper we suggest a method for direct evaluation of the average and reliability 
error exponent over the ensemble of LDPC error- correcting codes of given rate and 
connectivity. An analytical solution has been obtained, for both Gallager and MN 
codes, using methods of statistical physics, which is in perfect agreement with known 
results in the limit k,j{, t) oo (with R finite). The results for MN and Gallager codes 
become identical in this limit as both become random codes. 

Average error exponent results obtained by our method for codes of finite {k,j) 
values cannot be obtained using traditional approaches used in the information theory 
community. As expected, they seem to be upper bounded by the k,j^oo curves, but 
suggest a profoundly different behaviour for Gallager and MN LDPC codes. Average 
error exponent results for Gallager codes show a gradually improved performance as the 
parameters {k,j) increase, until they finally coincide with the k^j^oo error exponent 
result. The results for MN codes becomes identical to the k,j oo error exponent 
result for all > 2 or t > 2. To some extent, this is in agreement with previous 
results obtained for the critical flip rate of MN codes jS IH IH] and is a result of the 
close-to-random codebook they generate. 

An interesting feature of the present study is the similarity of our equations to 
those obtained in [6^ in spite of the different approaches used. An important advantage 
offered by the current approach is a potential extension to select high performance codes 
to obtain reliability exponent Yahxes for LDPC codes of finite connectivity; obtaining such 
solutions remains a difficult task and is currently under study. 
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